Deconvolution apparatus and method using a local signal-to-noise ratio

ABSTRACT

A method for a deconvolution of a digital input image (I(x i )) having a plurality of input voxels (x i ), in particular a digital input image obtained from a medical observation device, such as a microscope or endoscope and/or using fluorescence, includes computing a local signal-to-noise ratio (SNR(x i )) within an input region (R(x i )) of the digital input image, the input region consisting of a subset of the plurality of input voxels of the digital input image and surrounding the current input voxel. A noise component (β(SNR)) is computed from the local signal-to-noise ratio, the noise component representing image noise ( ([h*f](x i ), n(x i )) in the deconvolution. The noise component is limited to a predetermined minimum noise value (β min ) for a local signal-to-noise ratio above a predetermined upper SNR threshold value (SNR max ) and is limited to a predetermined maximum noise value (β max ) for a local signal-to-noise ratio below a predetermined lower SNR threshold value (SNR min ).

CROSS-REFERENCE TO PRIOR APPLICATIONS

This application is a U.S. National Phase Application under 35 U.S.C. §371 of International Application No. PCT/EP2019/051863, filed on Jan.25, 2019, and claims benefit to European Patent Application No. EP18194617.9, filed on Sep. 14, 2018. The International Application waspublished in English on Mar. 19, 2020 as WO 2020/052814 under PCTArticle 21(2).

FIELD

The invention relates to a method and an apparatus for digital imagerestoration, in particular digital input images from microscopes orendoscopes.

BACKGROUND

In the area of image processing, deconvolution is known to enhance theresolution and reduce the noise in a digital input image. The underlyingassumption of deconvolution is that the observed digital input imageI(x_(i)) results from the true image {circumflex over (f)}(x_(i)) asfollows:I(x _(i))=

([h*{circumflex over (f)}](x _(i)))+n(x _(i)),wherein h(x_(i)) designates a linear transfer function, e.g. from therecording system, such as a point-spread function, n(x_(i)) designatesadditive image noise,

designates Poisson noise depending on the signal strength, x_(i)designates an input voxel, and * denotes a convolution. In thedeconvolution, i.e. the computation of an estimate of the true image{circumflex over (f)}(x_(i)), the image noise is represented by a noisecomponent.

The digital input image may be two or higher dimensional. Thus, putgenerally, the input image may be n-dimensional, n being an integerlarger than one, n≥2.

At each input voxel x_(i), the input image is represented by at leastone image parameter I, a digital value, which in particular may berepresentative of light intensity. A voxel may be a pixel if the digitalinput image is two-dimensional. More generally, a voxel may be a ndimensional data structure if the input image is n dimensional. Then-dimensionality may result from e.g. a spatial—three-dimensional—imageand/or a two-dimensional image having more than one color channel. Forexample, the input image may comprise at each voxel a plurality ofintensity values at different spectral bands, such as any one orcombination of an R, G, B value and/or other colors in the visible andnon-visible-light range, UV, IR, NIR values, or a single grey scaleintensity value.

Here and in the following, an input voxel marks a location x_(i) in theinput image, where x_(i) may be represented by the n-tuple x_(i)={x₁, .. . , x_(n)} of n local coordinates for a n dimensional image. Forexample, in a two-dimensional RGB image, n=5 holds, because there aretwo spatial coordinates x, y and three color channels of the image, R,G, B: x_(i)={x₁, . . . , x₅}={x, y, R, G, B}. Alternatively, an RGBinput image may be considered as three independent two-dimensional greyscale images R(x, y), G(x, y), B(x, y), where R, G, B denote e.g. thelight intensity I in the respective color band. The input image may alsobe three-dimensional. For a three-dimensional grey scale image, n=3,x_(i)={x₁, x₂, x₃}={x, y, z}.

For the deconvolution of a noisy input image, knowledge about thesignal-to-noise ratio is needed. In practice, the image noise and thesignal-to-noise ratio must be estimated. Thus, the deconvolution allowsonly computing an approximation f(x_(i)) of the true input image{circumflex over (f)}(x_(i)).

SUMMARY

In an embodiment, the present invention provides a method for adeconvolution of a digital input image (I(x_(i))) having a plurality ofinput voxels (x_(i)), in particular a digital input image obtained froma medical observation device, such as a microscope or endoscope and/orusing fluorescence. A local signal-to-noise ratio (SNR(x_(i))) within aninput region (R(x_(i))) of the digital input image is computed. Theinput region consists of a subset of the plurality of input voxels ofthe digital input image and surrounding the current input voxel. A noisecomponent (β(SNR)) is computed from the local signal-to-noise ratio, thenoise component representing image noise (

([h*f](x_(i)), n(x_(i))) in the deconvolution. The noise component islimited to a predetermined minimum noise value (β_(min)) for a localsignal-to-noise ratio above a predetermined upper SNR threshold value(SNR_(max)) and is limited to a predetermined maximum noise value(β_(max)) for a local signal-to-noise ratio below a predetermined lowerSNR threshold value (SNR_(min)).

BRIEF DESCRIPTION OF THE DRAWINGS

Embodiments of the present invention will be described in even greaterdetail below based on the exemplary figures. The present invention isnot limited to the exemplary embodiments. All features described and/orillustrated herein can be used alone or combined in differentcombinations in embodiments of the present invention. The features andadvantages of various embodiments of the present invention will becomeapparent by reading the following detailed description with reference tothe attached drawings which illustrate the following:

FIG. 1 shows a schematic representation of a medical observationapparatus according to an embodiment of the invention;

FIG. 2 shows an example of an input region;

FIG. 3 shows a schematic flowchart of a method according to anembodiment of the invention;

FIG. 4 shows an input image before deconvolution;

FIG. 5 shows an output image computed from the input image of FIG. 3using a known deconvolution method;

FIG. 6 shows an output image computed from the input image of FIG. 3using the deconvolution according to an embodiment of the invention.

DETAILED DESCRIPTION

Standard deconvolution algorithms fail to deliver optimum results inspecific applications, such as microscopy or endoscopy, in which theimage characteristics in particular of fluorescence images differ fromthose encountered in other applications such as portrait, landscape,entertainment or artistic imagery. Moreover, an elaborated estimationand computation of the signal-to-noise ratio may quickly become tooburdensome, in particular in view of the increased spatial resolution ofmodern images, to be performed in real time, especially if deconvolutionhas to keep up with a high frame rate, such as 60 Hz.

Embodiments of the present invention provide an apparatus and methodwith an improved deconvolution performance in particular in case offluorescence, endoscope and/or microscope images. Further, embodimentsof the present invention provide an apparatus and method, which are ableto deconvolve a high-resolution digital input image in real time.

According to an embodiment of the invention the improved deconvolutionis provided by a method, being related to a method for a deconvolutionof a digital input image having a plurality of input voxels, inparticular a digital input image obtained from a microscope or endoscopeand/or using fluorescence, the method comprising the steps of:

computing a local signal-to-noise ratio within an input region of thedigital input image, the input region consisting of a subset of theplurality of input voxels of the digital input image and surrounding thecurrent input voxel;

computing a noise component from the local signal-to-noise ratio, thenoise component representing image noise in the deconvolution;

wherein the noise component is limited to a predetermined minimum noisevalue for a local signal-to-noise ratio above a predetermined upper SNRthreshold value and is limited to a predetermined maximum noise valuefor a local signal-to-noise ratio below a predetermined lower SNRthreshold value.

Further, the improved deconvolution is provided according to anembodiment of the invention by an image processor for a medicalobservation apparatus, the image processor comprising:

a storage or memory section for storing a digital input image comprisinga plurality of input voxels and

a deconvolution engine for computing a deconvolved output image from theplurality of input voxels;

wherein the deconvolution engine comprises a noise component, whichdepends on a local signal-to-noise ratio at an input voxel, the localsignal-to-noise ratio being computed in an input region consisting of asubset of the plurality of input voxels of the digital input image;

wherein the image processor contains a predetermined upper SNR thresholdvalue and a predetermined lower SNR threshold value; and

wherein the noise component is limited to a predetermined minimum noisevalue for a local signal-to-noise ratio above the predetermined upperSNR threshold value and to a predetermined maximum value for a localsignal-to-noise ratio below the predetermined lower signal-to-noisethreshold value.

The solution according to embodiments of the invention results in animproved quality i.e. reduced noise and increased sharpness, of thedeconvolved output image by using a local signal-to-noise ratio, whichis computed for the input region only. At the same time, the inventivemethod and apparatus according to embodiments realize a computationaladvantage by limiting the computation of the additive-noise component tothe local input region surrounding the input voxel and by introducingupper and lower SNR threshold values.

One reason of the improved quality of the output images obtained by theinventive method and apparatus according to embodiments may be thatparticularly in images of biological matter, such as obtained bymicroscopes and endoscopes and/or by using fluorescence, thecharacteristics of the noise may differ across the entire digital inputimage. For example, the noise characteristics in a region of the image,which is in focus and shows an illuminated or fluorescing part of a cellor of tissue may not be the same as for a region containing mostly blackor white background.

Moreover, it has been observed that for a large local signal-to-noiseratio above the predetermined upper SNR threshold, the image is goodenough, and no further computation of the regularization parameter needsto be performed. Thus, for a signal-to-noise ratio above thepredetermined upper SNR threshold, the noise component can be set to apredetermined minimum noise value.

Further, it has been observed that for very a small localsignal-to-noise ratio, SNR<1 or even SNR<<1, allowing the noisecomponent to increase beyond a predetermined maximum value does not leadto better results in the deconvolved output image. Thus, the noisecomponent is limited to the predetermined maximum SNR value if thesignal-to-noise falls below a predetermined lower SNR threshold.

The solution according to embodiments of the invention may be furtherimproved by adding one or more of the following features, which can becombined independently from one another.

For example, the image processor may be a software device, a hardwaredevice, or a combination of both a hardware device and a softwaredevice. The image processor may comprise at least one of a CPU, an arrayprocessor, a GPU, an ASIC, a FPGA and/or a FPU. The image processor mayin particular be adapted to perform array and/or parallel processingassembly. The image processor may comprise one or more software modulesthat, in operation, alter the physical structure of a processor by e.g.altering transistor states.

The signal-to-noise ratio SNR(x_(i)) at a voxel x_(i) of the input imagemay be computed in the input region R(x_(i)) as the ratio of signallevel or strength S(x_(i)) to noise level N(x_(i)):

${{SNR}\left( x_{i} \right)} = {\frac{S\left( x_{i} \right)}{N\left( x_{i} \right)}.}$

In another embodiment, the signal-to-noise ratio may be computed as:

${{SNR}\left( x_{i} \right)} = {10\mspace{11mu}\log_{10}{\frac{S\left( x_{i} \right)}{N\left( x_{i} \right)}.}}$

Once the local signal-to-noise ratio has been computed for all inputvoxels x_(i) of the digital input image, the resulting n-dimensionalarray SNR(x_(i)) may be blurred using a linear filter such as a low-passor a Gaussian filter to ensure smooth transitions between the voxels.

Preferably, only the pixels in a predetermined input region surroundingthe current voxel are used to compute the local signal-to-noise ratio atthat particular current voxel. In particular, the noise level N(x_(i))and the signal level S(x_(i)) may each be computed in the (same) inputregion only.

According to one embodiment, the signal level or strength S(x_(i)) at aninput voxel x_(i) may be estimated by convolving the input imageI(x_(i)) with a blur kernel k_(b)(x_(i)) in the input region R(x_(i))for the local SNR estimation and taking the maximum value of thisconvolution:

${{S\left( x_{i} \right)} = {\max\limits_{\Omega\;{R{(x_{i})}}}\left\lbrack {\left\lbrack {I*k_{b}} \right\rbrack\left( x_{i} \right)} \right\rbrack}},$where Ω is the region where the convolution is carried out.

According to one embodiment, the input region may contain between 200and 1000 input voxels, more preferably between 300 and 700 input voxelsand most preferably around 500 input voxels. An input region of thissize is sufficiently large to give reliable statistical estimates and atthe same time is not too large to extend to regions in an image, whichhave a different noise characteristics.

The blur kernel may be a linear filter, such as a Gaussian filter, e.g.having the following form:

${{k_{b}\left( x_{i} \right)} = {\frac{1}{\left( {2\pi} \right)^{n/2}{\prod\limits_{i = 1}^{n}\;\sigma_{i}}}e^{- {\sum\limits_{i = 1}^{n}{\frac{1}{2}{(\frac{x_{i}}{\sigma_{i}})}^{2}}}}}},$wherein 0.5<σ_(i)<1.0, in particular σ_(i)≈0.75. The linear filter mayhave e.g. a dimension of 3×3 or 5×5. On a GPU, directly computing such atwo-dimensional blur kernel is faster than separate one-dimensionalpasses using e.g. a 1×3 or 1×5 kernel.

According to another embodiment, the noise level N(x_(i)) at an inputvoxel x_(i) may be estimated in an advantageous embodiment by computinga variance based on the digital input image in the input region asfollows:N(x _(i))=√{square root over (var[I′(x _(i))])},

Herein, I′(x_(i)) is a derivative value or gradient of the input imageat location x_(i) in any of the n dimensions of the input image and ofany order. The derivative value I′ may for example be obtained byapplying a linear discrete operator

to the input image I in the input region R:I′(x _(i))=

I(x _(i)).

In particular,

may be an edge detection filter, such as a Sobel or Laplace operator, ora first, second or higher order gradient operator such as

$\frac{\partial}{\partial x},\frac{\partial}{\partial y},\frac{\partial^{2}}{{\partial x}{\partial y}},\frac{\partial^{2}}{\partial y^{2}},\mspace{14mu}{{or}\mspace{14mu}\frac{\partial^{2}}{\partial x^{2}}},$or any linear combination of such linear operators. For example, theapplication of the second-order gradient operator

${\mathbb{D}}\hat{=}\frac{\partial^{2}}{{\partial x}{\partial y}}$in the discrete may be represented in the two-dimensional case as:

I _(ab) =I _((a+1)(b+1)) −I _((a+1)b) −I _(a(b+1)) +I _(ab)at the input voxel x_(i)={a, b}, Ω∈R(x_(i)), and I_(ab)=I(a, b), where Ωis the considered volume or region.

In one embodiment, the variance may be computed as:var[I(x _(i))]=

[(I′(x _(i))−

[I′(x _(i))])²],where

is a norm of an n-dimensional array, such as a mean or expected valueoperator, for example

${{\mathbb{E}}(z)} = {\frac{1}{Q}{\sum\limits_{q = 1}^{Q}\; z_{q}}}$in the one-dimensional case, where z_(q) is a one-dimensional discretearray. As is explained further below, this form is particularly suitedfor computing the noise level in a computationally highly efficientmanner using a summed-area table.

For computing the (local) signal-to-noise ratio SNR(x_(i)), and forblurring the signal-to-noise ratio, the image processor may comprise asignal-to-noise computation engine, which may be a hardware or asoftware component, or a combination of both, of the image processor.The signal-to-noise computation engine may comprise in one embodiment aseparate hardware component, such as at least one CPU, GPU, FPU and/ordesignated IC, e.g. an ASIC, or any combination thereof, and/or aseparate software module, such as one or more subroutines. Thesignal-to-noise computation engine may in particular be configured tocompute also at least one of the signal level and the noise level.

As stated above, at least one of the signal level and the noise level atan input voxel x_(i), preferably both, may be computed in the inputregion R(x_(i)) only. The input region thereby consists of only a subsetof the input voxels of the digital input image.

The input region preferably is a contiguous region of input voxels. Theinput region preferably surrounds the current input pixel or voxel, forwhich the signal-to-noise ratio SNR(x_(i)) is computed. Thus, thesignal-to-noise ratio is a local value, which varies from input voxel toinput voxel, as different input regions, each of preferably the sameshape and size, are used for each input voxel.

Preferably, the input region, in which the input voxels are consideredfor computation of the signal and/or noise level, has at least two-foldsymmetry. The input region may deviate from an n-dimensional cuboidshape by a plurality of smaller n-dimensional cuboid regions, each ofwhich is contiguous. For example, the input region may differ in thetwo-dimensional case from a rectangle or a square shape by one or morerectangles that themselves are not part of the input region.

In particular, the input region may be a n-dimensional Gaussian,circular, Mexican hat or star-shaped contiguous arrangement, or adiscrete approximation thereof, of input voxels, in particular ofrectangles comprising a plurality of input voxels.

The input region may have a shape that results from superimposing aplurality of n-dimensional cuboid regions of input pixels. In thetwo-dimensional case, the input region may be constituted by asuperposition of rectangles, which each may comprise at least 2×2 or 4×4voxels.

According to one embodiment, a larger weight may be assigned to inputvoxels in the input region that are closer to the current input voxel atwhich the noise level is computed than to input voxels in the inputregion that are spaced further apart from the current input voxel. Sucha weighting leads to a smooth fading out of the input region towards itsborders and thus to a smoother distribution of the local noise levelsacross the input voxels.

The computation of the local noise level using

[(I′(x_(i))−

[I′(x_(i))])²] within the input region is preferably done using asummed-area table. This allows determining the local signal-to-noiseratio in a computationally highly efficient manner. A description ofsummed-area tables is found in Crow, Franklin C.: “Summed-area tablesfor texture mapping”, ACM SIGGRAPH computer graphics, (18) 3, 1984. Ann-dimensional version of a summed-area table is described in Tapia,Ernesto: “A note on the computation of high-dimensional integralimages”, Pattern Recognition Letters, 32 (2), 2011.

For each parameter that needs to be summed up over the input region, aseparate summed-area table may be generated, preferably before the noiselevel is computed.

For example, when computing the variance as

[(I′(x_(i))−

[I′(x_(i))])²], a summed-area table may be generated for

[I′_(i)(x_(i))] and, once this summed-area table has been computed, aseparate summed-area table may be generated for

[(I′(x_(i))−

[I′(x_(i))])²].

The at least one summed-area table may be part of the image processor,e.g. by being stored at least temporally in a storage section of theimage processor and/or by an instruction allocating memory of the imageprocessor for storing a summed-area table.

Defining the noise level in terms of a variance, which in turn iscomputed using a summed-area table, leads to a computationally highlyefficient deconvolution algorithm, which decreases considerably theadditional burden that the computation of a plurality of localsignal-to-noise ratios instead of a single global signal-to-noise ratiocreates.

Each of the summed-area tables is a digital array of values having thesame dimensionality as the input image.

Whereas it is described above that the noise level may be estimated bycomputing the variance, i.e. the second statistical moment, of thederivative value in the input region, it is also possible to use ahigher, n^(th), statistical moment for computing the noise level. Ingeneral, n^(th) the statistical moment:N(x _(i))=

[(I′(x _(i))−

[I′(x _(i))])^(n)]may also be computed using summed-area tables.

Using summed-area tables, the size and shape of the input region may beimplicitly defined by the selection of coordinates within thesummed-area table relative to the current input voxel, which coordinatesdefine the blocks, i.e. n-dimensional cuboids, of summed values that areused for computing the statistical moment. By combining the blocks ofthe summed-area tables to construe the input region, input voxelslocated where the blocks overlap automatically obtain a larger weight.It is in such a case preferred that the overlap of the blocks is locatedat the center of the input region. Thus, this region may automaticallybe assigned a larger weight.

In case of the image processor according to the invention, the imageprocessor may comprise a summed-area-table generator, thesummed-area-table generator being configured to compute a summed-areatable of the region and to compute the noise level using the summed-areatable. The summed-area table generator may be a hardware device, asoftware device, or a combination of both a hardware and a softwaredevice.

To further increase computational efficiency, at least one summed-areatable of the input region may be generated by a parallel computation ofprefix sums. In order to perform this step quickly, the image processor,or in particular the summed-area table generator preferably comprises atleast one prefix-sum generator, which may be a hardware device, asoftware device or a combination of the two. The hardware device maycomprise an array processor, such as a GPU, or a parallel-processingassembly comprising e.g. a plurality of GPUs, CPUs and/or FPUs. Thesoftware device may comprise an instruction set for the hardware device.In operation, the hardware set may change its structure due to theinstruction set of the software device to carry out the parallelcomputation of the prefix sums. The prefix-sum generator is configuredto compute a prefix sum in a parallel manner.

The deconvolution may be done using different methods.

For example, a Wiener deconvolution may be carried out using a Wienerfilter. In this embodiment, the output image f(x_(i)) can be estimatedfrom the observed image, i.e. the digital input image as follows:

${{f\left( x_{i} \right)} = {\mathcal{F}^{- 1}\left\lbrack \frac{{\mathcal{F}\left\lbrack {h\left( x_{i} \right)} \right\rbrack}^{*}{\mathcal{F}\left\lbrack {I\left( x_{i} \right)} \right\rbrack}}{{{\mathcal{F}\left\lbrack {h\left( x_{i} \right)} \right\rbrack}}^{2} + \beta} \right\rbrack}},$where

denotes the n-dimensional discrete Fourier transform of the quantity inbrackets and

[h(x_(i))]* is the complex conjugate of the discrete Fourier transformof h(x_(i)). The term:

$\mathcal{W} = \frac{{F\left\lbrack {h\left( x_{i} \right)} \right\rbrack}^{*}}{\left| {F\left\lbrack {h\left( x_{i} \right)} \right\rbrack} \middle| {}_{2}{+ \beta} \right.}$is the Wiener filter. The term β is the noise component of the Wienerfilter, which in one preferred embodiment is chosen to depend on thelocal signal-to-noise ratio SNR as described above: β=β(SNR(x_(i))).

In another embodiment, a maximum a posteriori deconvolution, such as aLucy-Richardson deconvolution, may be used. The Lucy-Richardsondeconvolution is e.g. described in Fish D. A.; Brinicombe A. M.; Pike E.R.; Walker J. G.: “Blind deconvolution by means of the Richardson-Lucyalgorithm”, Journal of the Optical Society of America A, 12 (1), p.58-65 (1995). The Lucy-Richardson deconvolution may in particular usescaled-gradient projection for increased computation performance.

In the Lucy-Richardson deconvolution, the following iterativecomputation of the output image is performed:

${{f^{({k + 1})}\left( x_{i} \right)} = {\frac{f^{(k)}\left( x_{i} \right)}{1 + {\beta{V\left( x_{i} \right)}}}{h^{T}\left( x_{i} \right)}*\left( \frac{l\left( x_{i} \right)}{{\left\lbrack {h*f^{(k)}} \right\rbrack\left( x_{i} \right)} + {b\left( x_{i} \right)}} \right)}},$with f^((k))(x_(i)) as the output image, i.e. the estimate of the trueimage, in the k-th iteration. Convolution is denoted with * andh^(T)(x_(i))=h(−x_(i)) is the flipped PSF (point spread function) of theimaging system or the lens. The function V(x_(i)) is the derivative ofthe regularization function, R_(reg)(x_(i)):

${V\left( x_{i} \right)} = {\frac{\delta R_{reg}}{\delta{f^{(k)}\left( x_{i} \right)}}.}$

The regularization function R_(reg)(x_(i)) may be at least one of thetotal variation ∥∇f(x_(i))∥², i.e. the norm of the gradient, a Tikhonovregularization ∥f(x_(i))∥² and Good's roughness

${\frac{\nabla{f\left( x_{i} \right)}}{f\left( x_{i} \right)}}^{2},$or any preferably linear combination of such functions.

The term b(x_(i)) is a background parameter. Preferably, the backgroundparameter is chosen to be dependent on the local signal-to-noise ratio:b(x_(i)):=b(SNR(x_(i))). According to a further embodiment, thebackground parameter may be a linear function of the localsignal-to-noise-ratio.

In an embodiment that leads to a particularly sharp deconvolved image,the background parameter attains a minimum value if the localsignal-to-noise ratio is above the predetermined lower SNR thresholdvalue. Additionally or alternatively, the background parameter mayattain a maximum value if the local signal-to-noise ratio is below thepredetermined upper SNR threshold value.

The background parameter may include a scaling factor that depends onthe noise level as e.g. computed above for the entire input image.

The parameter β corresponds to the noise component of the deconvolutionand is chosen in one embodiment of the invention to depend on the localsignal-to-noise ratio. β:=β(SNR(x_(i))). The local signal-to-noise-ratiomay be computed as described above.

Irrespective of the particular deconvolution method applied to thedigital input image, the noise component β(SNR) is computed as acontinuous function of the local signal-to-noise ratio. In particular,the noise component may be a monotonical, preferably strictlymonotonical function decreasing between a predetermined maximum noisevalue, β_(max), at the predetermined lower SNR threshold value,SNR_(min), β(SNR_(min))=β_(max) and a predetermined minimum noise value,β_(min), at a predetermined maximum SNR value, SNR_(max).

In particular, the noise component β may comprise at least one of alinear function, a polynomial function and a trigonometric functiondepending on

$\frac{1}{SNR}.$

According to one specific embodiment, which leads to good deconvolutionresults, β_(max) may be chosen to be between 1 and 10% of the dynamicrange, in particular about 5%. In the latter case, β_(max) may be about0.05. The value for SNR_(min) may be zero, i.e. SNR_(min)=0.

In particular for digital input images recorded by a microscope or anendoscope, the predetermined maximum SNR value, SNR_(max), at which theregularization parameter is truncated, is larger than 2, preferablysmaller than 5 and most preferred about 3.

Good results for deconvolution have been obtained if the gradient of thenoise component,

$\frac{d{\beta\left( {SNR} \right)}}{dSNR},$is smaller at SNR=0 and/or at the predetermined maximum SNR valueSNR_(max) than in the middle region between those two values, i.e. atabout

$\frac{{SNR_{\max}} - {SNR}_{\min}}{2}.$

According to a further preferred embodiment, which leads to very gooddeconvolution results, the noise component β is computed using:

${{\beta\left( {SN{R\left( x_{i} \right)}} \right)} = {\max\left\lbrack {{SNR}_{\min^{\prime}},{\frac{\beta_{\max}}{2^{\prime}}\left\{ {1 - {\arctan\left( \frac{\left( {{2\frac{SN{R\left( x_{i} \right)}}{SNR_{\max}}} - 1} \right)\pi}{2} \right)}} \right\}}} \right\rbrack}},$where SNR_(max) is the predetermined maximum SNR value

The improved deconvolution is also provided by a medical observationdevice, such as a microscope or an endoscope, comprising an imageprocessor in at least one of the above embodiments.

Finally, another embodiment of the invention also relates to anon-transitory computer readable medium storing a program causing acomputer to execute the image processing method in any one of theabove-mentioned embodiments.

In the following, a practical implementation of the invention isdescribed with reference to the drawings using one or more embodiments.The description serves as an example only and in no way restricts theinvention and its features to what is shown in the drawings andexplained below. Rather, the combination of features shown in theembodiments may be altered as described above. For example, one or moreof the features of an embodiment can be omitted if its or theirtechnical effect is not needed for a particular implementation.Likewise, one or more of the above-described features may be added ifthe technical effect of that particular feature or features isadvantageous for a particular application.

In the drawings, the same reference numerals are used for elements thatcorrespond to each other with respect to at least one of structure andfunction.

First, the structure of an image processor 1 according to an embodimentof the invention is described with reference to FIG. 1. The imageprocessor 1 may comprise one or more hardware components 2, such as amemory 4 and an integrated circuit 6. The integrated circuit 6 maycomprise at least one CPU, GPU, FPU, ASIC FPGA and/or any combinationthereof. The image processor 1 may be part of a computer 8, such as apersonal computer, laptop or tablet. Instead of, or in addition to theone or more hardware components 2, the image processor 1 may comprisesoftware, which is configured to run on the one or more hardwarecomponents 2.

As shown exemplarily in FIG. 1, the image processor 1 may be part of amedical observation device 10, such as a microscope 12, in particular aconfocal scanning microscope, or an endoscope.

Using an optical recording device, such as a camera 14 having a lens 16,the medical observation device 10 records a digital input image I(x_(i))which comprises input voxels x_(i). The digital input image I(x_(i)) isrepresented by a n-dimensional array of digital values, which represente.g. intensity in a color channel or, equivalently, a spectral band. Forexample, a digital input image I(x_(i)) may consist of one or more colorchannels R, G, B in the case of an RGB camera.

The camera 14 may be a RGB-camera, a monochrome camera, a multispectralor a hyperspectral camera or an assembly comprising a combination ofsuch cameras, of which the images may be combined. The camera 14 may bein particular a camera, which is sensitive to the spectral band orspectral bands, in which a fluorophore emits fluorescence.

The digital input image I(x_(i)) may be two- or more dimensional,wherein both spatial and color dimensions may be considered in thedimensionality. Although the camera 14 may record a time-series 18 ofdigital input images I(x_(i)), a single digital input image I(x_(i)) isconsidered in the following. This single digital input image maynonetheless have been computed from more than one original input images,as is e.g. the case with images obtained by z-stacking or withHDR-images. The digital input image I(x_(i)) may also be a stereoscopicimage.

In the case of an endoscope, a lens 16 of the microscope 12 may bereplaced by fiber optics. Otherwise, the endoscope and the microscopemay be considered as being identical for the purpose of the invention.

The image processor 1 may comprise several components, of which each maybe configured as a hardware device, as a software device, or as acombination of both a hardware and a software device, like the imageprocessor 1 itself.

For example, the image processor 1 may comprise an input interface 20for inputting the digital input image I(x_(i)). The input interface 20may provide a wired or wireless connection 22 to the camera 14. Thedigital input image I(x_(i)) is sent via the connection 22 from thecamera 14 to the image processor 1. The input interface 20 may provideone or more standard connections for receiving the digital input image,such as a HDMI, DVI, Bluetooth, TCP/IP, or RGB connection as well as anyother type of connection suited for the task of transmitting videoand/or image data.

The image processor 1 may further comprise an output interface 24 foroutputting a deconvolved digital output image f(x_(i)) to a storagemedium and/or one or more display devices 26, such as an eyepiece 28, amonitor 30 and/or virtual-reality glasses 32. Any of the one or moredisplay devices 26 may be stereoscopic or holographic.

The one or more display devices 26 may be connected wired and/orwireless to the output interface 24. The output interface 24 may employone or more data transmission standards for sending the deconvolveddigital output image f(x_(i)) to the one or more digital displaydevices, such as HDMI, DVI, Bluetooth, TCP/IP, or RGB, or any otherstandard that allows sending digital image data.

The image processor 1 is configured to perform a deconvolution of thedigital input image I(x_(i)) to compute the deconvolved digital outputimage f(x_(i)).

In the embodiment shown, the image processor 1 comprises a deconvolutionengine 34, which is configured to perform the deconvolution of thedigital input image I(x_(i)). The deconvolution may be a Wienerdeconvolution or a maximum a posteriori deconvolution, in particular aLucy-Richardson deconvolution.

The image processor 1 may further comprise a summed-area-table generator36, which is configured to compute at least one summed-area table 38 ofthe digital input image I(x_(i)) using a parallel prefix sum algorithm.The summed-area-table generator 36 may, as the image processor 1, be ahardware device, a software device, or a combination of both.

In the embodiment shown, the summed-area table generator 36 isconfigured to compute a summed-area table of both

[

I(x_(i))] and of

[(

I(x_(i))−

[

I(x_(i))])^(n)] for the entire input image. Here,

is a linear derivative operator, such as a Sobel or a Laplace filter ora gradient of any order. In the embodiment shown,

$\hat{=}\frac{\partial^{2}}{{\partial x}{\partial y}}$and n=2,

denotes the mean of an array.

From the summed-area tables, a background parameter b₀ of the overallimage can be computed quickly as follows:b ₀(x _(i))=0.2√{square root over (

_(Ω)[(

I(x _(i))−

_(Ω)[

I(x _(i))])²])}∀(x _(i) ∈I(x _(i))).

The expression

_(Ω) is used to denote that the operation is carried out only in theregion Ω which is equal to or a subset of the input region R(x_(i)). Forcomputing the local signal-to-noise ratio at a voxel x_(i) of the inputimage, an input region R(x_(i)) is defined. The input region surroundsthe voxel x_(i) at which the signal-to-noise ratio SNR(x_(i)) is to becomputed at a later step. The input region is composed of rectangularblocks in case of a two-dimensional input image and of cuboids in caseof a n-dimensional input image. The input region R(x_(i)) containsbetween 300 and 700 voxels, preferably about 500 voxels.

In the exemplary embodiment, the input image is two-dimensional and theinput region is constituted of rectangles. Preferably, a plurality ofrectangles overlap at the center of the input region R(x_(i)), inparticular at the current input voxel x_(i). Such an overlap leadsautomatically to a larger weight of the center region when a summed-areatable is used. The increased weight at the center of the input regioncreates a smooth tapering out at the edges of the input region and thusto better results.

An example of an input region R(x_(i)) centered on the current voxelx_(i) resulting from the superposition of three rectangles 40, 42, 44 isshown in FIG. 2. The input region may have at least two-fold symmetry.Of course, other shapes of input regions may also be used. It ispreferred that the input region R(x_(i)) approximates a circular region.For increased computational efficiency in using the summed-area table,the input region should contain as few rectangles and/or as largerectangles as possible. By construing the input region from asuperposition of rectangles, a weighing takes place automatically. InFIG. 2, the weighing factors are indicated in the different rectangularparts of the input region. Parts closer to the center of the inputregion, where the current voxel is located have a higher weight thanparts farther away from the center. In the example shown, the weighingfactors are 3, 2 and 1, respectively.

The noise level N(x_(i)) at an input voxel x_(i) is then computed byusing only the voxels in the input region R(x_(i)). The shape of theinput region R(x_(i)) is determined implicitly by defining therectangles 40, 42, 44 of which the sum is extracted from the summed-areatables. The noise level is computed as follows:N(x _(i))=0.2√{square root over (

_(Ω)[(

I(x _(i))−

_(Ω)[

I(x _(i))])²])}∀(x _(i) ∈R(x _(i))).

Next, the signal level S(x_(i)) is computed by again using only thevoxels in the input region R(x_(i)) as follows:

${S\left( x_{i} \right)} = {\max\limits_{\Omega = {R{(x_{i})}}}\left\lbrack {\left\lbrack {I*k_{b}} \right\rbrack{\left( x_{i} \right).}} \right.}$

The expression k_(b)(x_(i)) denotes a two-dimensional Gaussian blurkernel:

${{k_{b}\left( x_{i} \right)} = {\frac{8}{9\pi}e^{{- \frac{8}{9}}\mspace{11mu}{\sum\limits_{i = 1}^{1}x_{i}^{2}}}}}.$

Of course, other blur kernels may also be used.

Once both the noise level N(x_(i)) and the signal level S(x_(i)) havebeen determined for each of the input voxels x_(i), the signal-to-noiseratio is determined for each input voxel x_(i) as:

${SN{R\left( x_{i} \right)}} = {\frac{s\left( x_{i} \right)}{N\left( x_{i} \right)}.}$

The resulting digital array SNR(x_(i)) has the same dimensionality asthe input image I(x_(i)). Preferably, the array SNR(x_(i)) is low-passfiltered or, more preferably, blurred with a 3×3 or a 5×5 Gaussianfilter with σ=0.8.

If it is chosen to use e.g. a Lucy-Richardson deconvolution, thefollowing iteration is carried out to compute the deconvolved outputimage f(x_(i)):

${{f^{({k + 1})}\left( x_{i} \right)} = {\frac{f^{(k)}\left( x_{i} \right)}{1 + {{\beta\left( {SN{R\left( x_{i} \right)}} \right)}{V\left( x_{i} \right)}}}{h^{T}\left( x_{i} \right)}*\left( \frac{J\left( x_{i} \right)}{{\left\lbrack {h*f^{(k)}} \right\rbrack\left( x_{i} \right)} + {b\left( x_{i} \right)}} \right)}},$until the following convergence criterion has been met:

$\frac{\sum\limits_{{over}\mspace{14mu}{all}\mspace{11mu}{elements}}{{{f^{(k)}\left( x_{i} \right)} - {f^{({k - 1})}\left( x_{i} \right)}}}}{\sum\limits_{{over}\mspace{14mu}{all}\mspace{11mu}{elements}}\left( {{f^{(k)}\left( x_{i} \right)} + {f^{({k - 1})}\left( x_{i} \right)}} \right)} < {2 \times 10^{- 5}}$

Here, f^((k))(x_(i)) represents the deconvolved output image at thek^(th) iteration. The noise component β(SNR(x_(i))) is computed as

${{\beta\left( {SN{R\left( x_{i} \right)}} \right)} = {\max\left\lbrack {{SNR}_{\min^{\prime}},{\frac{\beta_{\max}}{2^{\prime}}\left\{ {1 - {\arctan\left( \frac{\left( {{2\frac{SN{R\left( x_{i} \right)}}{SNR_{\max}}} - 1} \right)\pi}{2} \right)}} \right\}}} \right\rbrack}},$wherein SNR_(max) represents a maximum SNR threshold value, which is setto 4 in the current embodiment.

The background parameter b(x_(i)) is chosen to be a linear function ofthe signal-to-noise ratio in the current embodiment. It is computed as

${b\left( x_{i} \right)} = {{\max\left\lbrack {0,\ {b_{0}\left( {1 - \frac{SN{R\left( x_{i} \right)}}{SNR_{\max}}} \right)}} \right\rbrack}.}$

Finally, the expression V(x_(i)) is computed as

${V\left( x_{i} \right)} = {{\frac{\delta}{\delta{f\left( x_{i} \right)}}{\frac{\nabla{f\left( x_{i} \right)}}{f\left( x_{i} \right)}}^{2}} = {{- 2}{\nabla\left( \frac{\nabla{f\left( x_{i} \right)}}{f\left( x_{i} \right)} \right)}}}$

If the convergence criterion has been met, f^((k+1))(x_(i)) representsthe digital deconvolved output image. The output image f^((k+1))(x_(i))is a digital array of values and has the same dimensionality as theinput image I(x_(i)).

The output image f^((k+1))(x_(i)) is then output to the output interface24 from where it may be passed on to any of the one or more displaydevices 26 for display and/or a storage device.

The deconvolution carried out by the image processor 1 may be summarizedas shown in FIG. 3.

At the start of the method, constants and functions used in thedeconvolution may be selected by a user, or automatically, e.g. aspreset values, by the image processor 1.

For example, the upper SNR threshold value, SNR_(max), the lower SNRthreshold value, SNR_(min), may be defined or altered.

The form of the noise component β(SNR(x_(i))) may be selected by a usere.g. from a list of available functions.

The specific form of the operators

and/or

and/or the order of the statistical moments for computing the noiselevel may be selected by a user.

The preferably linear filter used for blurring the array SNR(x_(i))and/or the blur kernel k_(b)(x_(i)) may be selected by a user e.g. froma list of available filters.

The particular deconvolution method may be selected by a user. Forexample, the user may select between a Wiener deconvolution and amaximum a posteriori deconvolution such as a Lucy-Richardsondeconvolution.

If a Wiener deconvolution is selected, the user may specify or select atransfer function, e.g. from a library, such as a point-spread functionappropriate for the current recording system. Alternatively, thetransfer or point-spread function may be stored in a memory or storagesection and be automatically determined from the set-up of themicroscope or endoscope, e.g. by automatically detecting the type oflens which is currently in use.

If a Lucy-Richardson deconvolution is used, the regularization functionR_(reg)(x_(i)) may be selected by a user e.g. from a list of availablefunctions, such as the total variation ∥∇f(x_(i))∥², i.e. the norm ofthe gradient, a Typhonov regularization ∥f(x_(i))∥² and Good's roughness

${\frac{\nabla{f\left( x_{i} \right)}}{f\left( x_{i} \right)}}^{2}$or any combination thereof. The function used in the backgroundparameter b(x_(i)) may be selected by a user, e.g. by selecting afunction from a list.

A user may define the shape and size of the input region R(x_(i)), e.g.by selecting from a list of available input regions or by specifying anindividual input region.

Further, any constant used in the deconvolution may be set and/oraltered by the user.

In a step 50, the digital input image I(x_(i)) is acquired. This stepmay be carried out by the camera 14 or by the image processor 1, whichmay e.g. compute the digital input image I(x_(i)) from one or moredigital images residing in the memory 4, acquired from the camera 14, ora combination of both.

In the next step 52, the at least one summed-area table 38 is computedfrom the input image I(x_(i)) preferably using a parallel algorithm suchas a prefix sum. Step 52 may be carried out at any time before thedeconvolution is computed. As indicated above, two summed-area tables 38are computed if the noise level is computed using a statistical momenthaving an order of at least 2.

Once the summed-area tables 38 are available, the local noise levelN(x_(i)) is computed for every voxel x_(i) in the input image I(x_(i))in step 54, using only the voxels in the input region R(x_(i)) aroundthe current voxel x_(i).

At any time before computing the deconvolution, the local signal levelS(x_(i)) is determined for every input voxel x_(i) of the input imageI(x_(i)) as described above, e.g. using a Gaussian blur kernel, andindicated at step 56. For the computation of the signal level S(x_(i))at any input voxel x_(i) of the input image I(x_(i)), only the voxelsx_(i) in the input region R(x_(i)) are used.

Once the signal level S(x_(i)) and the noise level N(x_(i)) areavailable, the signal-to-noise ratio is determined for every input voxelx_(i) of the input image at step 58.

The resulting array SNR(x_(i)) containing the signal-to-noise ratio foreach voxel x_(i) may be blurred using a low-pass and/or a Gaussianfilter. This is denoted by step 60.

At step 62, the deconvolution is carried out and the output imagef(x_(i)) is computed

The output image f(x_(i)) may undergo further post-processing as isindicated at step 64. For example, the output image may be assignedpseudo-colors and/or be merged with other images if the input imageI(x_(i)) is a fluorescence image in the emission band of a fluorescingfluorophore.

At step 66, the output image f(x_(i)) may be displayed on a displaydevice 26 and/or stored in a storage device, such as a disk or memorycard.

The results of a deconvolution of an input image I(x_(i)) becomeapparent from FIG. 4 to 6. In FIG. 4, input image I(x_(i)) is shown. InFIG. 5, an output image is shown as obtained by a conventionalLucy-Richardson deconvolution in which only the global signal-to-noiseratio for the entire image was used in the background parameter and thenoise component.

In FIG. 6, the result of using the deconvolution as described in thecontext of the exemplary embodiment of the invention is shown.

A comparison of FIG. 5 and FIG. 6 clearly shows that the deconvolutionof the embodiment according to the invention is capable of reducingnoise more efficiently and of rendering finer detail than theconventional deconvolution.

At the same time, defining the noise level in terms of a statisticalmoment and of summed-area tables computed with prefix sums allows aquick computation of the deconvolved output image f(x_(i)) although, inprinciple, the use of a local signal-to-noise ratio is computationallymore tedious than using a global signal-to-noise ratio.

While embodiments of the invention have been illustrated and describedin detail in the drawings and foregoing description, such illustrationand description are to be considered illustrative or exemplary and notrestrictive. It will be understood that changes and modifications may bemade by those of ordinary skill within the scope of the followingclaims. In particular, the present invention covers further embodimentswith any combination of features from different embodiments describedabove and below. Additionally, statements made herein characterizing theinvention refer to an embodiment of the invention and not necessarilyall embodiments.

The terms used in the claims should be construed to have the broadestreasonable interpretation consistent with the foregoing description. Forexample, the use of the article “a” or “the” in introducing an elementshould not be interpreted as being exclusive of a plurality of elements.Likewise, the recitation of “or” should be interpreted as beinginclusive, such that the recitation of “A or B” is not exclusive of “Aand B,” unless it is clear from the context or the foregoing descriptionthat only one of A and B is intended. Further, the recitation of “atleast one of A, B and C” should be interpreted as one or more of a groupof elements consisting of A, B and C, and should not be interpreted asrequiring at least one of each of the listed elements A, B and C,regardless of whether A, B and C are related as categories or otherwise.Moreover, the recitation of “A, B and/or C” or “at least one of A, B orC” should be interpreted as including any singular entity from thelisted elements, e.g., A, any subset from the listed elements, e.g., Aand B, or the entire list of elements A, B and C.

REFERENCE NUMERALS

-   1 image processor-   2 hardware component-   4 memory-   6 integrated circuit-   8 computer-   10 medical observation device-   12 microscope-   14 camera-   16 lens-   18 time-series of input images-   20 input interface-   22 connection between camera and image processor-   24 output interface-   26 display device-   28 eyepiece-   30 monitor-   32 virtual reality glasses-   34 deconvolution engine-   36 summed-area-table generator-   38 summed-area table-   40, 42, 44 rectangles constituting the input region-   50 step of acquiring the digital input image-   52 step of computing the summed-area table-   54 step of computing the noise level-   56 step of computing the signal level-   58 step of computing the signal-to-noise ratio-   60 step of blurring the signal-to-noise ratio array-   62 step of computing the output image by deconvolving the input    image-   64 post-processing of the output image-   66 step of displaying and/or storing the output image-   b₀, b(x_(i)) background parameter-   f(x_(i)) deconvolved digital output image-   f^((k+1))(x_(i)) (k+1)^(th) iteration in an iterative computation of    the output image-   {circumflex over (f)}(x_(i)) true image-   h(x_(i)) transfer function of recording system-   I(x_(i)) digital input image-   k_(b)(x_(i)) blur kernel-   N(x_(i)) noise level-   R(x_(i)) input region for determining at least one of the signal    level, the noise level and the signal-to-noise ratio-   R_(reg)(x_(i)) regularization function-   S(x_(i)) signal level-   SNR(x_(i)) signal-to-noise ratio-   SNR_(max) (predetermined) upper SNR threshold value-   SNR_(min) (predetermined) lower SNR threshold value-   V(x_(i)) functional derivative of regularization function-   x_(i) voxel-   derivative operator-   mean of an array-   a, b, c, n, x, y coordinates, variables-   β(x_(i)) noise component of the deconvolution-   β_(max) (predetermined) maximum noise value-   β_(min) (predetermined) minimum noise value-   Ω region of the input image, in which an operation is carried out

The invention claimed is:
 1. A method for a deconvolution of a digitalinput image (I(x_(i))) having a plurality of input voxels (x_(i)), inparticular a digital input image obtained from a medical observationdevice, such as a microscope or endoscope and/or using fluorescence, themethod comprising: computing a local signal-to-noise ratio (SNR(x_(i)))within an input region (R(x_(i))) of the digital input image, the inputregion consisting of a subset of the plurality of input voxels of thedigital input image and surrounding the current input voxel, andcomputing a noise component (β(SNR)) from the local signal-to-noiseratio, the noise component representing image noise (

([h*f](x_(i)), n(x_(i))) in the deconvolution, wherein the noisecomponent is limited to a predetermined minimum noise value (β_(min))for a local signal-to-noise ratio above a predetermined upper SNRthreshold value (SNR_(max)) and is limited to a predetermined maximumnoise value (β_(max)) for a local signal-to-noise ratio below apredetermined lower SNR threshold value (SNR_(min)).
 2. The methodaccording to claim 1, wherein the step of computing the localsignal-to-noise ratio (SNR) includes computing a local signal level(S(x_(i))) in the input region (R(x_(i))) and a local noise level(N(x_(i))) in the input region (R(x_(i))), and wherein the step ofcomputing the noise level includes computing at least one summed-areatable for at least one contiguous region of input voxels (x_(i)).
 3. Themethod according to claim 2, wherein the step of computing the at leastone summed-area table includes the step of computing the at least onesummed-area table in a parallel manner using a prefix sum.
 4. The methodaccording to claim 2, wherein computing the local noise level (N)includes computing a variance of the digital input image (I(x_(i)))using at least one summed-area table.
 5. The method according to claim4, wherein the variance is computed after applying a linear derivativeoperator (

) to the input image data (I(x_(i))) to obtain derivative image data(I′(x_(i))).
 6. The method according to claim 5, wherein the linearderivative operator (

) is a gradient operator or an edge-detection filter.
 7. The methodaccording to claim 1, wherein the step of computing the localsignal-to-noise ratio (SNR(x_(i))) includes the step of computing asignal level (S(x_(i))) at an input voxel (x_(i)) in the input region(R(x_(i))), and wherein the step of computing the signal level includesconvolving the input image data (I(x_(i))) with a blur kernel(k_(b)(x_(i))).
 8. The method according to claim 1, wherein the noisecomponent (β(SNR(x_(i)))) is computed as having a gradient$\left( \frac{{d\beta}({SNR})}{dSNR} \right),$ which is smaller at thepredetermined upper SNR threshold value (SNR_(max)) and/or at thepredetermined lower SNR threshold value (SNR_(min)) than between thepredetermined upper SNR threshold value (SNR_(max)) and thepredetermined lower SNR threshold value (SNR_(min)).
 9. The methodaccording to claim 8, wherein the step of computing the noise component(β(SNR)) includes the step of computing a trigonometric function. 10.The method according to claim 1, wherein the deconvolution is aLucy-Richardson deconvolution.
 11. A non-transitory computer readablemedium storing a program causing a computer to execute the imageprocessing method according to claim
 1. 12. A medical observationapparatus, such as a microscope or endoscope comprising an imageprocessor, the image processor being configured to carry out the methodaccording to claim
 1. 13. An image processor for a medical observationapparatus, such as a microscope or endoscope, the image processorcomprising: a memory configured to store a digital input image(I(x_(i))) comprising a plurality of input voxels (x_(i)), and adeconvolution engine configured to compute a deconvolved output image(f(x_(i))) from the plurality of input voxels, wherein the deconvolutionengine comprises a noise component (β(SNR(x_(i)))) which depends on alocal signal-to-noise ratio (SNR(x_(i))) at an input voxel, the localsignal-to-noise ratio being computed only in an input region (R(x_(i)))consisting of a subset of the plurality of input voxels of the digitalinput image, wherein the image processor contains a predetermined upperSNR threshold value (SNR_(max)) and a predetermined lowersignal-to-noise threshold value (SNR_(max)), and wherein the noisecomponent is limited to a predetermined minimum noise value (β_(min))for a local signal-to-noise ration above the predetermined upper SNRthreshold value (SNR_(max)) and to a predetermined maximum value(β_(max)) for a local signal-to-noise ratio below the predeterminedlower signal-to-noise threshold value (SNR_(min)).
 14. The imageprocessor according to claim 13, wherein the image processor comprises asummed-area-table generator, the summed-area-table generator beingconfigured to compute a summed-area table of the digital input image andto compute the local noise level using the summed-area table.
 15. Amicroscope or endoscope comprising the image processor according toclaim 13.